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In forced polymer translocation, the average translocation time, r, scales with respect to pore force, /, and 
polymer length, N, as r ~ f~ 1 N 13 . We demonstrate that an artifact in Metropolis Monte Carlo method 
resulting in breakage of the force scaling with large / may be responsible for some of the controversies between 
different computationally obtained results and also between computational and experimental results. Using 
Langevin dynamics simulations we show that the scaling exponent j3 < 1 + v is not universal, but depends on 
/. Moreover, we show that forced translocation can be described by a relatively simple force balance argument 
and P to arise solely from the initial polymer configuration. 



I. INTRODUCTION 

The force-driven transport of biopolymers through a nano- 
scale pore in a membrane is a ubiquitous process in biol- 
ogy. Despite the complex dynamics involved in the pro- 
cess, Monte Carlo (MC) has been almost the only computa- 
tional method used for modeling it, see the schematic Fig.Q] 
Only fairly recently more realistic dynamics have been ap- 
plied GISIBS. 

The classic theoretical treatment based 
on writing down the free energy for a system of two equilib- 
rium ensembles separated by a wall JauIL henceforth called 
Brownian translocation, has not found support from MC sim- 
ulations. The close-to-equilibrium dynamics would validate 
the derivation of the translocation dynamics using the Rouse 
relaxation time or diffusion based arguments. The validity of 
assumptions about the diffusion constant along the polymer 
chain was questioned already by the authors of [8, 9], who 
also noted that the characteristic translocation time both for 
phantom and self-avoiding polymer chains was bound to be 
greater than their characteristic time for relaxation to thermal 
equilibrium. The theoretical treatment of forced translocation 
can be said to have so far been almost solely guided by MC 
simulations dEl, despite the MC results contradicting the 
available experimental results. Hence, the theory has evolved 
independently of the experimental findings. 

In the attempt to determine the dynamical universality of 
translocation processes the scaling of the average transloca- 
tion time t with respect to the polymer chain length N as 
t ~ N" has been under intensive study. There is an abun- 
dance of research reporting different scaling exponents, f3. 
Some of the lately reported results do find some consistency 
on (3. For example, both the MC and Langevin dynamics 
(LD) simulations were reported to give (3 = 2v for short and 
(3 — l + v for long polymers in 2D, where v is the swellin g ex - 
ponent lfl2l[l3ll . Strictly universal scaling was claimed in llllll . 
where j3 — 1.5 was obtained both in 2D and 3D. 

In addition, several computational investigations have sug- 
gested nearly linear scaling ((3 = 1). On one hand, this close- 
to-linear behavior is prone to appear at the over-damped {i.e. 
Brownian) limit without self-avoiding effects, as shown by 
MC simulations in II 1411 . On the other hand, in 11511 (3 value 
close to unity could be inferred, but the polymer chains were 
found to be clearly out of equilibrium both on cis and trans 



sides in Brownian dynamics simulations, which was ignored 
in the interpretation of the results. Linear scaling has also 
been obtained from MC simulations in 01611 and from an LD 
simulation 1 17], albeit for moderate ranges of N. 

Not only do the different computational results and accom- 
panying theories contradict, but they also seem to have very 
little to offer in explaining the available experimental results. 
For example, the experimental study by Storm et al. lfl8ll re- 
ports (3 = 1.27, whereas the FFE framework ITolfTlll claims 
uniform scaling different from this at all forces and is thus ob- 
viously in trouble here. Since there logically exists a require- 
ment for the results obtained by using more realistic computa- 
tional methods to be validated not only against experimental 
results but also against results from MC model simulations, a 
demand is placed for identifying the limitations and possible 
errors of the MC based translocation models. 

Our recent study using molecular-dynamics-based simu- 
lations showed that (3 depends on the pore force and that 
the obtained values for f3 are in accord with the experi- 
ments |5]. As pointed out by Storm et al [18], the exper- 
imental pore force magnitudes are anticipated to be larger 
than those used in computer simulations, which would ex- 
plain the differing f3 values. However, larger force magni- 
tudes in MC have been claimed to produce the same scaling 
relations as smaller ones. We suggest that saturation of transi- 
tion probabilities for large forces in MC could be responsible 
for this discre pancy b etween simulated and experimental (3 

values SdMEllHlil]. 

Accordingly, we set out to study the forced polymer translo- 
cation with the pore force, /, as a control parameter. Three 
different methods will be used for modeling the process. The 
standard MC using Metropolis sampling will be compared 
with kinetic (or n-fold) Monte Carlo method (KMC) llM |25ll 
in ID. This comparison is made in order to characterize the 
forced translocation process. Differences between MC and 
KMC results would indicate subtleties in the effective inde- 
pendence of the events comprising the system l26ll . LD sim- 
ulations will be used as a reference in ID. With LD we also 
compute results in 2D and 3D to gain further evidence for our 
description of the forced polymer translocation. Finally, we 
determine the effect of the initial configuration on the scaling 
of the translocation time, r. In order to be free from potential 
artifacts due to spatial discretization, all results are produced 
by off-lattice simulations. 
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FIG. 1: Snapshots from 3D Langevin dynamics simulations, a) A 
straight initial configuration of a chain with N — 50 beads, b) A 
relaxed equilibrium initial configuration of a chain with N — 50 
beads. The wide solid line is the pore-to-end distance R pe . c) A 
chain with iV = 100 beads translocating with a force of //-Do = 
2.91, when s = 80 beads have translocated. 



II. BACKGROUND 

The investigation of the dynamics of forced translocation 
is subtly related to the assumption of close-to-equilibrium dy- 
namics. The classic treatments of the problem by Sung and 
Park la] and Muthukumar [7] assumed that a free energy, F s , 
could be written down for the two ensembles on either side of 
the wall, when s segments have been translocated: 



^ s =(l- 7 ) In [*(#-*)] - SfAz 



k R T 



(1) 



where / is the pore force, and 7 = 0.69 for a self-avoiding 
chain [7fl. By using the free energy one can write down 
the one-dimensional Langevin equation for translocation as a 
function of translocation coordinate, s, when entropic effects 



are assumed small: 

ms — —m^dFg/ds + r)(t) = — m£/Az + rj(t), (2) 

where m, £, rj(t), and Az are the bead mass, friction con- 
stant, random force, and pore length, respectively. From the 
Langevin equation the scaling of the translocation time with 
the pore force r ~ / _1 is straightforwardly obtained, where 
/ = fb/kgT and b are the dimensionless pore force and the 
Kuhn length, respectively. Although the form of the Langevin 
equation (Eq. (O) is identical in forced and unforced translo- 
cation, the connection between the free energy, F s , and the 
pore force, /, breaks down in forced translocation, where the 
pore force is a control parameter. The Langevin equation ap- 
proach yields r oc A^ 2 and t oc iV 1 as bounds for the translo- 
cation time scaling 10, El. Based on mere unhindered motion 
of polymers over potential difference A/i ~ / and the initial 
equilibrium configuration, Kantor and Kardar predicted that 
the translocation time should scale as [9] 



N 
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(3) 



By the computational method where polymer follows de- 
tailed molecular dynamics and the solvent coarse grained 
stochastic rotation dynamics, we have previously found that 
the translocating polymers are driven increasingly out of equi- 
librium on both sides of the pore under pore force magnitudes 
relevant in biological and experimental systems The dy- 
namics is then mainly determined by the force balance be- 
tween the drag force exerted on the mobile beads on the cis 
side and the constant pore force, fd = /. Additional contri- 
bution comes from the crowding of the polymer beads on the 
trans side due to relaxation towards equilibrium being slower 
than the rate at which new segments enter through the pore. 
On the cis side the rate at which polymer beads are set in mo- 
tion towards the pore was seen to be greater than the rate at 
which beads entered the pore, which indicates that the tension 
spreads along the polymer faster than what the polymer is able 
to relax towards equilibrium. By measuring the rate at which 
the tension spreads, an estimate could then be made for the 
scaling of the translocation time with the polymer length. We 
will present these observations in the Results section obtained 
from LD that, due to the absence of hydrodynamics, shows 
these characteristics even more clearly. 



III. POLYMER MODEL 

In the model system adjacent monomers are connected with 
anharmonic springs, described by the finitely extensible non- 
linear elastic (FENE) potential, 



Ufene — — ttR 2 m (l — - — )■ 
Z H 



(4) 



Here I is the length of an effective bond, which can vary be- 
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equilibrium distance at which the bond potential takes its min- 
imum value. Choosing l max as the unit length and R = 0.3, 
yields l m i n = 0.4, and Iq — 0.7. In the standard MC simula- 
tions k B T = 1 and the spring constant K = 40 k B T. KMC 
dynamics proved more susceptible to bond fluctuations than 
standard MC, so k B T = 0.1 and K/k B T = 400 were used 
(see next section). The FENE potential suffices in ID. 

For 2D and 3D, where LD was used, the Lennard-Jones 
(LJ) potential, 



U LJ = 4e 



'<t\ 12 /a 
. r / V r 

U LJ = 0, r > 2^ 6 a, 



r < 2- 1/6 a (5) 
(6) 



was used between all and FENE-potential, Eq. (0]), between 
adjacent beads. The parameter values were chosen to be 
e = 1.2, and a = 1.0 for the LJ, and io = 0, R = 1.5 = l ma x, 
K = 60/er 2 for the FENE potential. The used LJ poten- 
tial with no attractive part mimics good solvent condition for 
the polymers. Initial states are relaxed equilibrium configura- 
tions. 



IV. TRANSLOCATION MODELS 

In ID we study the translocation dynamics using MC, 
KMC, and LD methods. In order to link the ID dynamics with 
forced translocation we define the segment s as translocated, 
when it passes the original position of the first bead. Time 
units cannot be expected to be identical for the three methods, 
but the scaling laws given by MC and KMC must be identical 
to those obtained from LD, if the MC and KMC methods are 
to be taken as representative models for forced translocation. 

One-dimensional off-lattice Monte Carlo simulations are 
performed with the standard Metropolis acceptance test. A 
transition to a new state is attempted by picking at random a 
bead, computing the change in the particle's potential energy 
resulting from an attempted trial move by a distance 5 from 
its present position r, AU = U(r + S) — U(r). Transition 
is always accepted for AU < 0, and according to its Boltz- 
mann weight exp(— (3AU), if AU > 0. Time is incremented 
after each attempt, whether accepted or rejected, by a constant 
amount At. For a detailed pseudocode, see Appendix A. 

Unlike in the Metropolis MC method, in KMC algorithm 
the system is moved to another state at every attempt, regard- 
less of how improbable the transition is. Accordingly, KMC 
has been the choice for doing simulations at low tempera- 
tures, where transition probabilities are low and, accordingly, 
Metropolis MC prohibitively slow. On the other hand, the ad- 
ditional bookkeeping required in KMC makes it computation- 
ally slower than MC at higher temperatures. However, MC 
and KMC have profound differences in other aspects than just 
computational efficiency. 

In KMC the probability of the move is reflected on the 
(stochastically) estimated time, At, for the transition to take 
place, see Appendix B. Estimation of the elapsed time in- 
volves computing the system's all transfer probabilities, or 
rates, pj, from which a cumulative function Ri — J2] =1 Pj 



is calculated. Ri is a class including events, whose proba- 
bilities pj G {pY im iPi iax ]- F rom m e possible transitions j 
the event to take place is picked at random, so that R4-1 < 
uR{t) < Ri, where u 6 (0, 1] is a random number and 
R(t) = Ri and I the number of transition classes. The time 
elapsed between the previous and current event is estimated as 
At = —1/R(t) lnu', where u' ^ u is a random number. It is 
noteworthy that the stochastically determined time increment 
At is inversely proportional to the the total probability, i.e. the 
rate of change, R(t), of the system evolving in time, here the 
instantaneous polymer configuration. 

Due to finite distance for the trial move, 5, the Metropo- 
lis MC inevitably eliminates very improbable moves, which 
in the case of our ID polymer simulations are those stretch- 
ing or compressing polymer bonds far from their equilibrium 
lengths, Iq. Hence, as seen in our simulations, bond length 
fluctuations are larger in KMC than in MC. To maintain sta- 
bility we used a lower temperature in our KMC simulations. A 
mere scaling of the bond potential magnitude by reducing the 
temperature does not change the dynamic universality class of 
the system. 

Our LD algorithm was implemented according to |27j] . The 
LD translocation algorithm was used also in two and three di- 
mensions. In the algorithm, the pore force, /, is a free control 
parameter and not derived from free energy, as was done in 
Eq. (|2). In the initial states the bond lengths are equal to re- 
laxation distances and the first bead is in the middle of the 
pore. The pore diameter is 1.26 and length 36, where 6 = 1 
is the Kuhn length of the modeled polymer. We have used 
k B T = 1, £ = 0.73, and m = 16 parameters in the Langevin 
implementation. This yields Dq = k B T/£m w 0.086 for the 
one-particle self-diffusion constant in ID. 



V. RESULTS 

Our motivation for evaluating MC based methods in the 
context of forced translocation comes on one hand from dis- 
crepancies between scaling relations r ~ N 13 obtained from 
different MC simulations and on the other hand discrepancy 
between all MC simulations and experiments. We also look 
into the observed separate scaling regimes for short and long 
polymers. We base our evaluation of the methods mainly 
on the scaling of the translocation time with pore force and 
with polymer length. Since we take LD as the reference for 
the physical translocation system when hydrodynamic inter- 
actions are ignored, we first check these scaling relations us- 
ing LD. 

Fig. [2] a) shows the scaling of the translocation time with 
respect to the dimensionless pore force r ~ f a obtained from 
LD simulations in 2D and 3D. The thermal energy k B T was 
kept fixed in all simulations. We obtained a = —0.990 ± 
0.01 in 2D and a = -0.978 ± 0.02 in 3D. So, the scaling in 
these dimensions is the expected r ~ / _1 in accordance with 
Eq. Q, confirming that entropic effects are weak. 

Fig. [2] shows the scalings of the translocation time with re- 
spect to chain length, r ~ N 13 , for a frictional and friction- 
less pore in 2D in b) and in 3D in c). The scaling exponents 
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(3, identical for longer chains with frictional and frictionless 
pores, are shown in Table J] for different values of f/D . 
There is a fair agreement between our 2D scaling exponents 
and those reported previously [13]. Qualitatively, the change 
of (3 with N with the frictional pore agrees with the aforemen- 
tioned LD results where two scaling regimes were claimed. 
Accordingly, we obtain a smaller (3 for shorter chains. How- 
ever, a quite different behavior is seen in the case of a friction- 
less pore. At a large enough force there exists only one scaling 
identical to the common scaling obtained for long polymers 
with both pores. Hence the translocation times longer than 
what would be expected from this scaling seem to be due to 
friction in the confined region inside the pore. 

TABLE I: Langevin dynamics. Values of /3 obtained from Fig. [2] 
Here U = /Az/fcsT is the dimensionless pore potential. 

U f/D pjd = 2) /3(rf = 3) 
0.0036 0.014 1.52 ±0.08 1.40 ±0.04 
0.0075 0.03 - 1.40 ±0.02 

0.03 0.35 1.54 ±0.03 1.43 ± 0.02 
0.75 2.91 1.58 ±0.01 1.44 ±0.03 
3 11 1.63 ±0.04 1.47 ±0.05 
30 116 1.68 ±0.02 1.50 ±0.01 

The scaling exponents (3 obtained in 2D, while agreeing 
with the results in lfl3ll . seem to be in strong disagreement 
with the scaling (3 — 1.53 obtained by 2D MC simulations for 
an infinite potential with the pore length Az = 1. In addition, 
in another MC simulation in 2D (3 = 1.70 has been obtained 
for two different forces / = 1, 5 with Az = 3 lfl2tl . So, not 
only MC results contradict with LD results, but MC simula- 
tions mutually disagree. To check if this could be due to MC 
translocation model possibly belonging to a different dynamic 
universality class from the respective LD model we simulated 
the simplest possible case, i.e. ID forced translocation using 
KMC in addition to MC and LD. We obtain invariably (3 — 2 
for LD, MC, and KMC, in agreement with the previous MC 
result [9], so the computational MC model using Metropolis 
sampling seems to produce correct time dependence, although 
forced translocation can hardly be taken as a purely Poisso- 
nian process |[2g] . However, we do find a method-dependent 
artifact which explains the apparent discrepancy. 

Also in ID the translocation time obtained from LD sim- 
ulations is found to scale with the pore force as r ~ 
Fig. |3] a). MC and KMC models also give this scaling with 
pore force / < 1. However, for / > 1 the translocation time 
levels off to a constant value. Potentially, KMC might not be 
as sensitive to this artifact, as will be discussed in more detail 
later. However, KMC algorithm fails to run with so large a 
force, due to the above-explained higher probability of mov- 
ing a particle regardless of whether a bond breaks or not. In 
2D the breaking of the force scaling has been reported in MC 
forced translocation simulations in two 111211 and three IU6I1 di- 
mensions. Our simulations strongly imply that it is related 
to the saturation of the transition probabilities in MC at large 
force values inside the pore (see Appendix A, step 6). 

The discrepancy between the different (3 values obtained 



from MC simulations can be explained by this artifact related 
to Metropolis sampling. The infinite force for which (3 = 1.53 
was obtained |9J] takes the MC model to the plateau region 
(/ > 1), where / = oo is no different from / = 1 in that it 
gives identical velocity, v, to the translocating polymer, when 
the pore length, Az, is constant. However, by increasing Az 
while keeping the force acting per polymer bead, /, constant, 
does increase pore potential and thus the polymer velocity, 
v, which in turn, according to our simulations, increases (3 
(see Fig. [2] and Table [IJ. The pore potential, AU, over the 
pore is explicitly shown to be the primary control parame- 
ter for forced translocation in Table [III where translocation 
times are compared when either the force acting on polymer 
beads inside the pore or the pore length is changed in the LD 
model. In U Az = 1, while in |Qj] Az = 3 was used. Of 

TABLE II: Results from 3D Langevin dynamics for a polymer of 
length N = 100. Here AU ps fAz/k B T is the dimensionless 
potential over the pore. 

AU f Az r 
30 10 3 326 ±10 
15 10 1.5 642 ±20 
30 20 1.5 354 ±10 

the two force values / = 1, 5 in lfl2ll the latter is within the 
plateau regime and the first at least very close to it, which 
would explain the obtained almost identical (3 values. Due 
to / = oo being in effect just the constant force value in the 
plateau regime, AU in |§] is lower than what was used in lfl~2tl . 
which explains why this (3 value for the infinite force is lower 
than what was obtained for force values 1 and 5. The satu- 
ration of the transition probabilities could also be the reason 
why in 3D MC simulation results presented in ITlll to support 
FFE predictions the value of (3 did not change, when the pore 
force was varied. Hence, it would seem that there is no reason 
to assume a universal (3 for all /. 

The prevailing discrepancy between simulations can be at 
least partly explained by the deficiency of the MC method at 
large force regime, which is the biologically and experimen- 
tally relevant force range. For example the Kuhn length of 
the single-stranded DNA (ssDNA) b — 1.6 nm, and a charge 
density 1.28 e/nm reported in [28] for free ssDNA in an ionic 
solvent would result in the effective charge of 2 e/b. Tak- 
ing into account that the pore is prone to screen the charges 
even more than the solvent ll29ll the effective charge on a 
nucleotide traversing the pore was estimated to be w 0.1 e. 
In the a-hemolysin pore of length 5.2 nm there are approx- 
imately 13 nucleotides of length 0.4 nm l30ll . yielding total 
effective charge q* = 1.3 e for the polymer segment inside 
the pore. Typical experimental pore potential is V « 120 
mV H d |M |H. At T = 300K, the dimensionless ra- 
tio U = q*V '/ksT w 6.03 is obtained, which corresponds to 
the dimensionless value U = fAz/ksT. Accordingly, the 
corresponding dimensionless force / = fb/ksT has a value 
of ~ 2 in our MC simulations, which is already well in the 
plateau regime, see Fig.[3]a). 

We checked the effect of the pore friction in ID, where dy- 
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FIG. 2: a) The scaling of the translocation time r with respect to the external force /. Results from 2D (A) and 3D (o) Langevin dynamics 
simulations yield the same / ~ iV^ 1 scaling. Since the simulations were done with the only difference in dimension, we note that in higher 
dimension the translocation process is faster, b) 2D and c) 3D Langevin dynamics simulation results for the scaling of the translocation time r 
with respect to chain length N with various forces. Open symbols are with frictional pore and full symbols are with a ballistic pore. All plots 
yield different values for the scaling exponent (3 which are reported in Table|l]along with their corresponding forces. The finite size effect due 
to the frictional is identical to ID case see text). With / = 10 only chains with N > 200 follow the asymptotic scaling behavior. If the pore 
is non-aqueous (i.e. ballistic) the finite size effect applies only for chains smaller than N = 30, and is amplified by reducing the force. 



namics is no more constrained inside than outside the pore, 
making comparison of the cases with frictional and friction- 
less pores straightforward. The same characteristics as in 
higher dimensions were seen. With a frictional pore (3 = 2 
is obtained only for TV > 30. From ID LD simulations 
with a frictional pore this maximum length, l max , was seen 
to increase with pore force, /. Increasing / increases the 
average velocity of the polymer, (v), and so increases fric- 
tional contribution. Since the force acts only on the part of 
the polymer inside the pore, at large enough force frictional 
term damps movement only in the direction from the cis to 
trans side thus increasing the translocation time proportion- 
ately more for short than long polymers. At a very weak force 
this effect is not perceived due to friction affecting movement 
in both directions (i.e. also from trans to cis), as seen also 
in dimensions 2 and 3 from the topmost curves in Figs. |2]b) 
and c). This is by definition a finite-size effect and bound to 
affect the obtained (3 more for shorter polymers, where the 
portion in the pore constitutes a larger part. The decrease of 
(3 due to pore friction is understandable, since in the asymp- 
totic limit of completely friction-dominated dynamics linear 
scaling, (3=1, with polymer length is to be expected. 

We measured the center-of-mass velocity, (v* m ), in the 
z-direction, which is the direction perpendicular to the wall 
from cis to trans, with a frictional and frictionless pore, see 
Fig. [3]b). Comparing the ID curves obtained from LD and 
MC in Figs. [3] a) and b) it is clear that (v* m ) is not simply 
reciprocal of r, which indicates that one must be cautious in 
applying equilibrium concepts and observables in characteriz- 
ing translocation. Based on this finding, using (zo)/t, where 
zq is the z coordinate of the last polymer bead in the initial 
equilibrated configuration, is not equivalent to using (v* m ) for 
translocation velocity. In ID it is particularly clear, that the 
first definition simply sets the translocation velocity v ~ 1 /r, 
a valid definition for translocation velocity, but one that also 
easily leads to drawing wrong conclusions about translocation 
dynamics. Comparing Fig.[2]a) and Fig. 0b), the same charac- 
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FIG. 3: a) Comparison of the scaling of translocation time r (in 
arbitrary units) with pore force, /, in ID obtained from MC (A), 
KMC ( x ) and LD (•) simulations. Polymers are of constant length 
N = 50. All computational methods yield the same force scaling 
r ~ / _1 for / < 1. The MC method shows an artifact in the force 
regime / ^> 1. KMC fails to perform in this regime (see text), b) 
The average center-of-mass velocity in the z-direction, («f m ), as a 
function of the dimensionless pore force /, for N = 20 (+, □) and 
N = 50 (x , o) in ID, for N = 100 (•) in 2D, and for N = 100 
(A) in 3D. The data for (□) and (o) was obtained using a frictionless 
pore in LD simulations. All other results are from LD simulations 
using a frictional pore. Linear scaling is plotted as a solid line to 
guide the eye. At the bottom, results from ID MC simulations (A) 
with iV = 50 are shown (in arbitrary units). See text for details. 



teristics is seen in dimensions 2 and 3. From Fig. 0b) it is seen 
first of all that the artifact in MC at large force magnitudes 
manifests itself in (v* m ) measured in the ID system. Outside 
of that, the (v£ m (/)) curves for ID MC and LD are qualita- 
tively fairly similar. The additional regime at very weak force 
magnitudes in LD was checked to be related with the form of 
the FENE potential, showing distinctively in ID where it is 
the only potential between polymer beads. 

To make sure that we are seeing finite-size effects for short 
polymers with a frictional pore instead of two scaling regimes 
for short and long polymers, we also check the dependence of 
(v*„, (/)} on polymer length, N. In accord with the findings 
in ||9|] the center-of-mass motion of a polymer translocating 
through a pore closely conforms to unimpeded motion, where 
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FIG. 4: a) The average center-of-mass motion along the transloca- 
tion coordinate, {vl m ), as a function of the chain length N in ID, 
with Langevin dynamics. Two different pore forces / = 0.1 (+, □) 
and 20 ( x , o) have been used. Both frictional (+, X) and frictionless 
(□, o) pores were used. MC results in ID with / = 0.1 (A) and 20 
(■) are also shown (in arbitrary units). The solid and dashed lines 
have a slope of — 1. b) (vl m ) as a function of the chain length N. 
Results are from Langevin dynamics simulations in 2D (□, o) and 
3D (+, x), with two different pore forces / = 0.25, 10 for the lower 
and higher points, respectively. Both frictional (+, □) and friction- 
less (x, o) pores were used. The solid lines ~ l/N are plotted to 
guide the eye. 
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FIG. 5: A snapshot at s = 35 of a 3D LD simulation with a polymer 
of length N = 100. The dimensionless pore force used here is / = 
0.01. The wide solid line represents the pore-to-end distance R pe as 
the number of mobile beads, s m , reaches the index of the observed 
bead, n = 52. 



( v cm) ~ 1/-^ with both the frictional and frictionless pore, 
see Fig. [4] It is evident that pore friction does not affect the 
behavior of (v* m ) that is a measure of the transfer of the whole 
polymer. Hence, the differences in Figs. |2]b) and c) have to 
come solely from the frictional contribution to the finite por- 
tion of the chain in the pore. The effect is seen clearly in r that 
is a local observable quantifying the transfer of the segments 
over the pore. 

The above observation that (uj m (/)) & l/r(f) indicated 
that the transfer of the part of the polymer crucial for charac- 
terization of forced translocation, i.e. the part in the vicinity 
of the pore does not remain even close to thermal equilibrium. 
Any meaningful description of the translocation process must 
then succeed to capture the mechanics of this part. Given that 
any diffusive, or close-to-equilibrium description is invalid, 
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FIG. 6: (Color online) a) 3D Langevin dynamics. Averaged squared 
distances of beads numbered 30, 50, 70„ 100, 200, and 400 from the 
pore as a function of the number of translocated beads s for polymers 
of length N = 100. The force / = 0.25. b) The number of mobile 
beads, s m vs the number of translocated beads, s, both normalized 
to the polymer length, N. 




0.25 0.5 0.75 



s/N 

FIG. 7: (Color online) 3D LD. Waiting time averages t(s) as a func- 
tion of the translocated segment s with f /Do = 2.91. Profiles with 
filled (N = 50: ■, N = 100: ; N = 400: A) and empty symbols 
(N — 50: □, N = 100: o) represent runs with an equilibrium initial 
configuration and with a straight initial configuration, respectively. 
Note that the latter profiles yield a clear linear behavior. 



the remaining minimum condition is force balance. Indeed, 
it turns out that the observed scaling behavior /3 < 1 + v 
in 2D and 3D can be explained by the force balance of the 
drag force exerted on the mobile beads on the cis side and 
the constant pore force, fd — f ■ Just as in we measured 
the squared distance, R^ e (n) (see Figs. Q]b) and|5]l, of the 
polymer bead n from the pore on the cis side during translo- 
cation, see Fig.|6]a), from which the segments toward the free 
end of the polymer are seen to remain immobile until they are 
pulled toward the pore. Fig. |6]b) shows that the number of 
mobile beads increases linearly as a function of translocated 
beads, s m = ks. Up to lengths ofJV« 100, k — N~ x , be- 
yond which it gradually levels off to a constant value greater 
than unity. The mobile beads comprise a moving segment, 
which already for moderate pore force exhibit no folding in- 
dicating that the relaxation of the beads in the moving seg- 
ment is far slower than the rate at which they are translocated, 
see Figs. [T]c) and [5] fd is exerted on the mobile beads, so 
fd ~ s m (v), where (v) is the average velocity of the mobile 
beads. When the whole chain has translocated, fd ~ N m (v), 
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where N m — kN . The beads are set in motion from their 
equilibrium positions, so the distance d of the last bead to be 
translocated scales as d ~ N v . The average translocation 
time then scales as r - (d)/(v) ~ kN 1+ " ~ N 1+V ~ x . For 
example for / /Dq = 2.9, we obtain \ ~ 0.3 from the data 
displayed in Fig. |6]b), which would give /3 w 1.3. The crowd- 
ing of polymer beads close to the pore on the trans side (see 
Fig- [He) increases (3 from this value. The crowding increases 
with force, which has been shown in dH]. From Fig. |6]b) it 
can be seen that k saturates when / is increased, which also 
increases (3 at large /. 

In order to confirm that the presence of v in the scaling ex- 
ponent (3 = 1 + v comes solely from the initial configuration 
and not from any factor due to close-to-equilibrium motion 
and Rouse relaxation time, we check the scaling and the aver- 
age waiting time, t(s), for the segment s to translocate for an 
initially extended configuration, Fig.Q]a). We obtain (3 = 2 
for this initial configuration. The average waiting time profiles 
of the initially extended and equilibrated configurations are 
shown in Fig. [7] For the initially extended configuration we 
obtain linear dependence t(s) ~ s giving exactly the obtained 
scaling with (3 = 2. Hence, the upper limit for the scaling 
exponent (3 results from the waiting times of the translocating 
beads, determined solely on their initial positions. In other 
words, already with a moderate pore force translocation veloc- 
ity completely dominates the process so that diffusive motion 
has no effect and the the only prevailing condition governing 
translocation dynamics is the force balance described above. 

VI. SUMMARY 

We have studied forced translocation without hydrodynam- 
ics in different dimensions using Monte Carlo (MC), kinetic 



Monte Carlo (KMC), and Langevin dynamics (LD) meth- 
ods. We have shown that forced translocation model using 
MC with basic Metropolis sampling gives the correct time 
dependence at moderate pore potentials but presents an arti- 
fact at large pore potentials. This artifact seems to explain 
the prevailing controversy between different MC results. It 
also seems to account for results claiming universal scaling 
of the translocation time with the polymer length, r ~ N@, 
independent of the pore force magnitude and hence the pre- 
vailing discrepancy between computational and experimental 
findings on forced translocation. 



Using LD in 2D and 3D we have shown that the scal- 
ing exponent (3 increases with the pore force, /. We obtain 
(3 = 1 + v as the high force limit. By measurements of poly- 
mers' center-of-mass velocities we have shown that descrip- 
tion of forced translocation with concepts related to close-to- 
equilibrium concepts is not well founded. We have given a 
description of the forced translocation based on simple force- 
balance condition, where the drag force exerted on the part of 
the polymer on the cis side changes with the number of poly- 
mer beads in motion. This description gives the above upper 
limit for (3 and also explains the increase of (3 with /. We 
have also shown that crowding of polymer beads is strong in 
the vicinity of the pore on the trans side. 

This work has been supported by the Academy of Finland 
(Project No. 127766). 
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VII. APPENDIX A : METROPOLIS MONTE CARLO 
ALGORITHM 

1. Choose an initial state, and set the time t = 0. 

2. Randomly choose a particle with label i, and calculate 
a trial position r- = fl + Srt. 

3. Calculate the energy change, AU = — /» ■ tJ% + 

U F ENE(n + Sri) - U FE NE(rl), resulting from this 
displacement. 

4. If AU < 0, the move is accepted; go to step 2. 

5. Get a uniform random number u 6 (0, 1]. 

6. If u < exp(— AU /UbT), accept the move. 

7. Update the time t = t + 1/N and go to step 2. 



VIII. APPENDIX B: KINETIC MONTE CARLO 
ALGORITHM 

The KMC algorithm for simulating the time evolution of a 
system where some processes occur with known average rates, 



or probabilities, pi = cxp(— AU), where AU is as defined in 
Appendix A, can be written as follows: 

1 . Choose an initial state, set the time t = 0, and form a 
list of all possible rates in the system pi. 

2. Calculate the cumulative function Ri = Y^j=iPj f° r 
i = 1, . . . , I, where / is the total number of transitions. 
Denote R~Ri. 

3. Get a uniform random number u 6 (0, 1]. 

4. Find the event to carry out by finding the i for which 

Ri-! <uR< Ri. 

5. Carry out event i. 

6. Recalculate all rates pi which may have changed due to 
the transition. If appropriate, remove or add new transi- 
tions i. Update N and the list of event rates accordingly. 

7. Get a new uniform random number v! G (0, 1]. 

8. Update the time with t = t + At, where At = - i2 fp 1 . 

9. Return to step 2. 

(note that the same average time scale can be obtained also 
using At = in step 8. However, including the random 
number describes better the stochastic nature of the process). 

This algorithm is known in different sources variously as 
the residence-time algorithm or the n-fold way or the Bortz- 
Kalos-Lebowitz (BKL) algorithm or just the kinetic Monte 
Carlo (KMC) algorithm. 



